﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace Nature.Toolkit.GKin
{
    /// <summary>
    /// 
    /// </summary>
    public abstract class EulerImplicitTPYSolver : TPYSoneSolver, IIdealGasMixtureTimeController
    {
        protected abstract void ComputeTP(ReadOnlyIdealGasMixture variable, TPYVector function);

        private IdealGasMixture m_basePoint;
        private double m_timeStep;
        private double m_time = 0.0;
        private double[] m_wyk;

        protected EulerImplicitTPYSolver(ModelObject model)
            : base(model)
        {
            m_basePoint = model.CreateIdealGasMixture();
            m_wyk = new double[model.Nsp];
        }

        protected override void Compute(ReadOnlyIdealGasMixture variable, TPYVector function)
        {
            var wk = variable.SpeciesMassFractionChemicalNetProductionRates;
            var yk = variable.SpeciesMassFractions;
            var yk0 = m_basePoint.SpeciesMassFractions;
            var speciesPart = function.Y;
            int nsp = wk.Length;
            double dt = m_timeStep; ;
            for (int isp = 0; isp < nsp; ++isp)
            {
                speciesPart[isp] = (yk[isp] - yk0[isp]) - dt * wk[isp];
            }
            this.ComputeTP(variable, function);
            m_time += dt;
        }

        public double MoveNext()
        {
            throw new NotImplementedException();
        }

        public virtual void Reset(IdealGasMixture mixture)
        {
            if (mixture.Model != Model) { throw new ArgumentException(); }
            m_basePoint.UpdateTPY(mixture);
            m_time = 0d;
            m_timeStep = 1.0e-11;
            m_basePoint.SpeciesMassFractionChemicalNetProductionRates.CopyTo(m_wyk);
            base.Solve(mixture);
            CalcTimeStep(
                m_basePoint.SpeciesMassFractionChemicalNetProductionRates.Data,
                mixture.SpeciesMassFractionChemicalNetProductionRates.Data);
        }

        void CalcTimeStep(double[] wk0, double[] wk1)
        {
            int nsp = wk0.Length;
            double dw = 0.0;
            for (int isp = 0; isp < nsp; ++isp)
            {
                double a = wk0[isp];
                double b = wk1[isp];
                double delta = Math.Abs(a - b) + double.Epsilon;
                double avg = 1.0e-10 * delta + Math.Max(Math.Abs(a), Math.Abs(b)) + double.Epsilon;
                double tmp = delta / avg;
                if (dw < tmp) { dw = tmp; }
            }
            m_timeStep *= 1.0e-2 / dw;
        }

        #region IIdealGasMixtureTimeController Members

        double IIdealGasMixtureTimeController.MoveNext() { return this.MoveNext(); }

        void IIdealGasMixtureTimeController.Reset(IdealGasMixture mixture) { this.Reset(mixture); }

        #endregion

        #region ICloneable Members

        object ICloneable.Clone()
        {
            throw new NotImplementedException();
        }

        #endregion
    }
}
